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We extend our previous work on applying CMB techniques to the mapping of gravitational- 
wave backgrounds to backgrounds which have non-GR polarisations. Our analysis and results 
are presented in the context of pulsar-timing array observations, but the overarching methods are 
general, and can be easily applied to LIGO or eLISA observations using appropriately modified 
response functions. Analytic expressions for the pulsar-timing response to gravitational waves with 
non-GR polarisation are given for each mode of a spin-weighted spherical-harmonic decomposition 
of the background, which permit the signal to be mapped across the sky to any desired resolution. 

We also derive the pulsar-timing overlap reduction functions for the various non-GR polarisations, 
finding analytic forms for anisotropic backgronnds with scalar-transverse (“breathing”) and vector- 
longitudinal polarisations, and a semi-analytic form for scalar-longitudinal backgrounds. Our results 
indicate that pulsar-timing observations will be completely insensitive to scalar-transverse mode 
anisotropies in the polarisation amplitude beyond dipole, and anisotropies in the power beyond 
qnadrupole. Analogously to our previous findings that pulsar-timing observations lack sensitivity to 
tensor-curl modes for a transverse-traceless tensor background, we also find insensitivity to vector- 
curl modes for a vector-longitudinal background. 
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I. INTRODUCTION 


A massive international effort is currently underway 
to observe gravitational waves across a wide range of 
frequencies. The second-generation of ground-based 
gravitational-wave interferometers are about to start col¬ 
lecting data, with Advanced LIGO [T] observation runs 
expected to begin before the end of 2015. The two Ad¬ 
vanced LIGO detectors will form part of a global net¬ 
work of kilometre-scale laser interferometers, with other 
instruments due to come online during the rest of this 
decade. These detectors will employ advanced technolo¬ 
gies to detect gravitational waves from stellar-mass com¬ 
pact binary systems emitting gravitational radiation in 
the kHz band EHS]. The European Space Agency re¬ 
cently selected a science theme based around a ~ 10® m 
arm-length space-based gravitational-wave interferome¬ 
ter (eLISA) for the L3 mission slot, due to launch in 
2034. Such a detector will observe gravitational waves 
in the millihertz band, which are generated by binaries 
involving the massive black holes that reside in the cen¬ 
tres of galaxies, with mass about one million times the 
mass of the Sun. These observations will permit tests 
of fundamental physics to exquisite precision, whilst also 
affording detailed demographic studies of massive black- 
hole populations [5]. 


Complementary to these experiments are ongoing 
efforts to characterise nanohertz gravitational waves 
through their perturbation to the arrival-times of radio 
signals from precisely timed ensembles of millisecond pul¬ 
sars spread throughout our galaxy [MH]- As a gravita¬ 
tional wave transits between the Earth and a pulsar, it 
induces a change in their proper separation, leading to a 
redshift in the arrival rate of the pulsar signals [IIHIl]- 
It is the exceptional stability of the integrated pulse pro¬ 
files of millisecond pulsars, and the resulting accuracy of 
the models for the pulse times of arrival (TOAs), that 
allow gravitational waves to be detected in this way. 

The differences between the modelled TOAs and the 
actual observed TOAs are known as the timing residuals. 
These residuals contain the influence of all unmodelled 
phenomena, such as additional receiver noise, interstel¬ 
lar medium effects, errors resulting from drifts in clock 
standards or ephemeris inaccuracies, and, most tanta- 
lisingly, gravitational radiation. The signature of grav¬ 
itational waves in these residuals may be deterministic 
or stochastic. The gravitational-wave sources expected 
to dominate the signal in the nanohertz frequency band 
are the early adiabatic inspirals of supermassive black- 
hole binary (SMBHB) systems [TSHTT] . Such systems 
are expected to form following the (suspected ubiqui¬ 
tous) mergers of massive galaxies during the hierarchical 
formation of structure. If there is a system which is par- 
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ticularly loud in gravitational-wave emission then this 
signal may be individually resolved and detected with 
pipelines dedicated to searches for the deterministic sig¬ 
nals of single sources [TB1120| . If, however, there are many 
sources which pile up in the frequency-domain beyond 
the ability of our techniques to separately resolve them, 
then the combined signal will form a stochastic back¬ 
ground of gravitational waves. Although there are other 
mechanisms which may contribute to a stochastic nHz 
gravitational-wave background (decay of cosmic-string 
networks [5IH21] or primordial remnants |25l 126)1. this 
incoherent superposition of signals from many SMBHB 
systems is expected to dominate the signal. 

Standard pipelines in use today employ cross¬ 
correlation techniques to search for stochastic back¬ 
grounds. The presence of a common background of grav¬ 
itational waves affecting the TOAs of all pulsars in an 
array (a so-called pulsar-timing array, PTA [27]) makes 
a cross-correlation search effective in leveraging the signal 
against uncorrelated noise processes. The concept of an 
overlap reduction function is common to stochastic back¬ 
ground searches for all types of gravitational-wave detec¬ 
tors, and describes the sky-averaged overlap of the an¬ 
tenna pattern functions of the two detectors whose data 
are being correlated |28|. In PTA analysis, the overlap re¬ 
duction function for a Gaussian, stationary, unpolarised, 
isotropic stochastic background composed of transverse- 
traceless (TT) gravitational-wave modes is a smoking- 
gun signature of the signal, known as the Hellings and 
Downs curve |5S|. It is a function of one variable: the 
angular separation between a pair of pulsars. 

For anisotropic distributions of gravitational-wave 
power on the sky, the overlap reduction function is 
no longer merely a function of the pulsars’ angular- 
separation. It will also depend on the positions of 
the pulsars on the sky relative to the distribution of 
gravitational-wave power, and thus will be a rich source 
of information in the precision-science-era of PTAs [501 
Furthermore, the overlap reduction function can be 
shown to vary when describing backgrounds where the 
graviton is permitted to have a small but non-zero mass 
[32j . The same is true when describing the overlap re¬ 
duction functions induced by gravitational-wave polari¬ 
sation states present in modified (metric) theories of grav¬ 
ity. In addition to the usual GR transverse-traceless ten¬ 
sor polarisation states, the beyond-GR polarisations con¬ 
sist of a scalar-transverse (“breathing”) state, a scalar- 
longitudinal state, and two vector-longitudinal states, 
each inducing correlation signatures which are markedly 
distinct from the Hellings and Downs curve [551 131 ] ■ 

In this paper we focus on the response of pulsar tim¬ 
ing observations to gravitational wave backgrounds with 
non-GR polarisation states. By decomposing a back¬ 
ground of given polarisation in terms of spin-weighted 
spherical harmonics, we are able to derive analytic ex¬ 
pressions for the detector response functions for each 
mode of each non-GR polarisation state as a function 
of the harmonic multipole. We discuss the implications 


of these results for mapping non-GR backgrounds to any 
desired angular resolution. We are also able to present 
analytic expressions for the overlap reduction functions 
of anisotropic scalar-transverse and vector-longitudinal 
backgrounds, whilst significant analytic headway is made 
for the corresponding function for scalar-longitudinal 
backgrounds. 

In Sec. [n] we introduce the concept of the measured 
signal in a gravitational-wave detector being a convolu¬ 
tion of the metric perturbations with the response ten¬ 
sor of the detector. We discuss the six distinct polar¬ 
isation states of gravitational waves which are permit¬ 
ted within a general metric theory of gravity by virtue 
of obeying Einstein’s Equivalence Principle. The ba¬ 
sis tensors for these polarisations are explicitly given. 
We also discuss the decomposition of the metric pertur¬ 
bations in terms of appropriate spin-weighted spherical 
harmonics. In [35] . the Fourier amplitudes of a plane- 
wave expansion of the metric perturbations for an arbi¬ 
trary transverse-traceless gravitational-wave background 
were decomposed in terms of a basis of spin-weight ±2 
spherical harmonics. In the case of scalar-transverse and 
scalar-longitudinal polarisations discussed in this paper, 
we decompose the Fourier amplitudes in terms of ordi¬ 
nary (spin-weight 0) spherical harmonics. For the vector- 
longitudinal polarisations, we decompose the Fourier am¬ 
plitudes in terms of spin-weight ±1 spherical harmonics. 
In Sec. [TTJ we also give expressions for the pulsar timing 
response functions, for either the polarisation or spin- 
weighted spherical harmonic expansion coefficients. The 
polarisation basis response functions for a pair of pulsars 
are given explicitly in the computational frame, where 
one pulsar lies along the z-axis and the other lies in the 
ccz-plane. These are needed for the overlap reduction 
functions calculations given in the following section. 

The overlap reduction functions for the different po¬ 
larisation states are studied in Sec. IIIII This function de¬ 
scribes the response of a pair of pulsars to a gravitational- 
wave background in a cross-correlation analysis, and is 
computed by integrating the overlap of the response of 
each pulsar to a particular gravitational-wave polarisa¬ 
tion over the entire sky. For a gravitational-wave back¬ 
ground with arbitrary angular structure, this sky inte¬ 
gral must be weighted by the gravitational-wave power 
at each sky location. We find an analytic expression 
for the overlap reduction function for a background with 
scalar-transverse (breathing) polarisation, and show that 
a PTA will lack sensitivity to angular structure beyond 
quadrupole in a cross-correlation analysis for this type of 
background. We also make significant analytic headway 
for the overlap reduction function of a scalar-longitudinal 
background, and find analytic forms for the limiting value 
in the case of co-directional and anti-directional pulsars. 
The overlap reduction function for a vector-longitudinal 
background with arbitrary angular structure is found an¬ 
alytically, with superficially perceived divergences in the 
overlap reduction function for co-directional pulsars re¬ 
solved by correctly incorporating the pulsar term in our 
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calculations. 

In Sec. |IV| we extend our previous work on mapping 
gravitational-wave backgrounds using CMB techniques 
[55] to non-GR polarisations. We derive analytic expres¬ 
sions for the response of a pulsar to each mode (cor¬ 
responding to a particular spin-weighted spherical har¬ 
monic) of the background, including the contribution 
from the pulsar term. In the process of doing these 
calculations, we find that the reason for the PTA in¬ 
sensitivity to angular structure beyond quadrupole in 
the gravitational-wave power of a scalar-transverse back¬ 
ground is due entirely to the corresponding lack of sensi¬ 
tivity of a single pulsar response to structure in the polar¬ 
isation amplitudes beyond dipole. We verify this analytic 
result with numerical map making and recovery. The pul¬ 
sar response to individual modes of a scalar-longitudinal 
and vector-longitudinal background are given analyti¬ 
cally, where in the latter case we find that PTAs com¬ 
pletely lack sensitivity to vector curl modes, analogous to 
our previous finding that PTAs lack sensitivity to tensor 
curl modes of a transverse-traceless background [35]. We 
discuss these findings further in Sec. jVj along with sug¬ 
gestions for future study and implications for the forth¬ 
coming analysis of real PTA data. 

Finally, we include several appendices (Apps. 
containing relevant information (e.g., definitions, identi¬ 
ties, recurrence relations) for spin-weighted and tensor 
spherical harmonics, Legendre polynomials, Bessel func¬ 
tions, etc., as well as providing technical details for the 
overlap reduction function and response function calcu¬ 
lations described in Secs. iniandllVl 


II. RESPONSE FUNCTIONS 


A. Detector response 


The response of a detector to a passing gravitational 
wave is given by the convolution of the metric perturba¬ 
tions hab{t,x) with the impulse response i?“*’(t,x) of the 
detector: 


r-(t) = f dr f d^yR°''’{T,y)hab{t-T,x-y). 
J —OO J 


( 1 ) 


If we write the metric perturbations as a superposition 
of plane waves 

habit, S)= r df [ ( 2 ) 

J — OO J S'^ 

then 

r{t) = r df [ d^O^ !?“''(/, k)habif, , (3) 

J — OO J 

where 


x[ dr (4) 


Further specification of the response function depends 
on the choice of gravitational-wave detector as well as on 
the basis tensors used to expand habif, k), as we explain 
below. 


B. Polarisation basis 

In standard GR, the Fourier components habif, h) are 
typically expanded in terms of the -I- and x polarisation 
basis tensors: 

habif, h) = h+if, fc)e+ (fc) -f /lx if, k)e^bi^) , (5) 

where 

eab(^) = ^Jb - Kk, 

^abi^) = 

and 0, (j) are the standard unit vectors tangent to the 
sphere: 

k = sin 9 cos (j)x + sin 9sm(j)y + cos 9 z , 

9 = cos9cos(j)x + cos9sm(j)y — sin9 z , (7) 

4> = — sin 4>x -\- cos (j) y. 

In this paper, we also consider modified metric theories 
of gravity, which admit four other types of polarisation: 
a scalar-transverse (or breathing) mode (R), a scalar- 
longitudinal mode (L), and two vector-longitudinal 
modes (A, Y). The polarisation basis tensors for these 
modes are: 


'abik) = Oak + ^ah, 

(8) 

•.abik) = V^kah , 

(9) 

'abih = ^ah + kaOb , 

(10) 

'abik^ — ^akb ka^b • 

(11) 


In terms of the polarisation tensors, the Fourier compo¬ 
nents habif, h) can be expanded generally as 

habif, k)=Y.hAif,k)eibik) (12) 

y4 

where A is some subset of {-f, x,i?,L, A, A}. The asso¬ 
ciated response function for a plane wave with frequency 
/, propagation direction k, and polarisation A is given 
by 

R^if,k) = R^\f,k)e^^ik), (13) 

and is related to the detector response r(/) via: 

/ OO P _ 

df d2U^5]i?^(/,fc)/iA(/,fc)e*""^‘. 

-OO 4 

(14) 

We will work with the polarisation basis response func¬ 
tions when calculating the various overlap reduction 
functions in Sec. imi 
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C. Spherical harmonic basis 


Eqs. (Bl) and (B9) in App. 


I 


Alternatively, we can expand the Fourier components 
habif, k) in terms of the appropriate spin-weighted spher¬ 
ical harmonics, as was done in [35]. A spin-weighted func¬ 
tion is a function of both position on the sphere, labelled 
k, and of a choice of an orthonormal basis, labelled l,m^ 
at points on the sphere. Under a rotation of the orthonor¬ 
mal basis, spin-weight functions transform in a particular 
way 


OO i 

hMk)=Y. E 


I —I rn— — l 


-\-(l 


Vc 

ilm) 



The above expressions for hab{f, k) can be written in 
compact form 


f{k, cos '0^ — sin0771, sini/jl + cosi/^m) = I, m) 


(15) 


^(/m) if')^(lm)ab (fc) (20) 

(Im) P 


where s is the spin-weight of the function. Any spin- 
weight s function can be expanded as a combination of 
spin-weighted spherical harmonics of the same weight, 
sYim{k). A spin-weight s spherical-harmonic can be re¬ 
lated to s derivatives of an ordinary spherical harmonic, 
as described in App. [^ 

For the standard GR tensor modes, if we define = 
Z“±7m“, we see that the combinations m'^m’^habif, k)aie 
spin-weight ±2 functions on the sphere. This allows the 
GR tensor modes to be expanded as combinations of 
spin-weight ±2 spherical harmonics, or equivalently in 
terms of the rank-2 gradient and curl sp heric al harmon¬ 
ics, y(?m)abCk), y(?m)abCk), defined by Fq. ^ in App.B 


OO L 

hab{f,k)=Y, E [^flm){f)y{?m)abCk) 


1—2 m— — l 


+^fl7n)if)y{lm)abik) ■ ( 16 ) 


if we take P to be a subset of {G, G, B, L, Vq, Vc}, and 
define 


y; 


\liabik) ^ ^yUk)e^i^Ck) (21) 


to unify the notation for the spherical harmonic basic 
tensors. (The factor of l/-\/2 is needed for the tensor 
spherical harmonics 1® satisfy orthonormality 

relations similar to Fqs. ( jBllj ) and (G8).) The associated 
response function for a given spherical harmonic mode is 

R(lm) if) = k)y^L)abik ), ( 22 ) 

and are related to the detector response r{t) via: 

/ OO 

(23) 

■°° (/m) P 


For the breathing and scalar-longitudinal modes, the 
functions m^rh^habif, k) are spin-weight 0 and so we can 
expand habif,k) in terms of ordinary (scalar) spherical 
harmonics: 

OO I 

K,if,k) = ^ E E ^fim)if)yUk)eZik), ( 17 ) 

1^0 m^-l 

- OO I 

habif,k) = -^Yl E ^\im)if)yi^ik)eabik), ( 18 ) 

y ^ 1=0 m=-l 


We will work with these response functions for the map¬ 


ping discussion in Sec. IV 


D. Pulsar timing response 

A gravitational wave transiting an Farth-pulsar line 
of sight creates a perturbation in the intervening metric. 
This causes a change in their proper separation, which is 
manifested as a redshift in the pulse frequency mHni: 


since the polarisation tensors e^^(fc) and e^f,(fc) are invari¬ 
ant under a rotation of 0, ip. For the vector-longitudinal 
modes, Th'^rh)}.hab{k) have spin-weight ±1 and so we 
can expand habif,k) in terms of spin-weight ±1 spher¬ 
ical harmonics or, equivalently, in terms of tensor fields 
y^m)abif^)’ y^m)abi^) Constructed from the rank-1 vec¬ 
tor spherical harmonics T(^j^(fc), T(^j^(fc) defined by 


z{t,k) 


Av{t) 

pq 




2l-hfc' 


Ahab{t,k) , 


(24) 


where k is the direction of propagation of the grav¬ 
itational wave, u is the direction to the pulsar, and 
Ahabit,k) is the difference between the metric pertur¬ 
bation at Farth, and at the pulsar some distance 

L from the Farth, (tp, Xp) = (t — L/c, x + Lu): 
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Ahab{t,k)= dfJ2hA{f,k)e^bik) 
J-OO 4 

/ oo _ 

d/ J2hA{f,k)e^,{k)< 

-OO A 


^i27Tf{t-k-xfc) _ ^i27rf{tp-k-Xp{c) 


i2-K f{t — k-x / c) 


1 — e 


— i27r/L(l+fc-n) / c 


(25) 

(26) 


For a gravitational wave background, which is a superposition of waves from all directions on the sky, the pulsar 
redshift integrated over k is given by 


J-oo Js^ ^ ^ 1 + fc • ^ 


Cab 


(fc) 


1 — e 


—i2T: f L{l-\-k-u) / c 




(27) 


Comparing the above expression with Eq. (14), we see that the detector response function R^{f,k) for a Doppler 
frequency measurement r(t) = z(t) is given by 


R^U,k) 


1 

2l + fc. 


u 


_ g-227r/L(l+fc-u)/c 


(28) 


For a timing residual measurement r{t) = f* dt'z{t'), the above response function R^{f, k) would need to be multiplied 
by a factor of ll{i2TTf). The response functions for individual spherical harmonic modes are similarly given by 


Rilm)if) 



1 

2 1 + fc • M 


{lm)ab 




2^ _ ^-i2TTfL(l+k-u)/c 


(29) 


E. Response functions for a pair of pulsars in the computational frame 


In the following section, we will calculate the correlated response of a pair of pulsars to a gravitational wave 
background. This calculation is most easily done in the so-called computational frame |30[ 1351136] , in which the two 
pulsars are in the directions 


ui = ( 0 , 0 , 1 ), 

U 2 = (sin ^,0, COS C) • 


(30) 


In addition, we can choose the origin of the computational frame to be at the solar-system barycentre (SSB), for 
which a detector (i.e., a radio telescope on Earth) has a: « (). In this frame the polarisation basis response functions 
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given in Eq. (28) simplify to 


Rtif,k) 

Rt{f,k) 

Rnf,k) 

R^{f,k) 

R?if,k) 

R2if,k) 

Ri{f,k) 

R^if,k) 

Riif,k) 

R2{f,k) 

Rlif,k) 

Rl{f,k) 


j(l - cos6») (l - , 


1 
2 

^\n a A a 2 sin^ C sin^ </< 

- (1 — sin C sin Cl cos (!) — cos Cl cos C)-:-^„ 

2 [ 1 + sin C sin 9 cos (j) + cos 9 cos C 


^ 27rz/_L2 (l+sin C sin 0 cos <ji+cos 0 cos ^)/c^ 


0 , 


1 / sin^ C cos 9 sin(2^) — sin(2C) sin 9 sin (p 


( I _ g-27ri/L2 (1+sin sin 9 cos ^+cos 9 cos (^)/c j 

/ ’ 

(1 - cos9) (l - e-^^^fLi{i+cose)/c^ ^ 

(1 - sin C sin 0 cos (^ - cos 0 cos C) (l - ® ® 


_ ^-27r2/Li(l+cos6»)/c^ 


_ g—27r2/L2 (1+sin sin 0 cos (/>+cos 0 cos ^)/c^ 

)■ 


1 cos^ 9 
yj2 1 + cos 9 

1 (sin^sin^cos^ + cos0cosC)^ 

•y/2 1 + sin sin 9 cos ip + cos 9 cos ( 

-cos6>sin6) / _ ^_ 2 ^i/Li(i+cose)/, 

1 + cos 9 V 

(sin Csin 9 cos (p + cos 9 cos C)(sin (p cos 9 cos cp — sin9 cos C) 

1 + sin sin 9 cos (p + cos 9 cos ( 

/-I —27r2/I/2 (1+sin sin 0 COS (?!)+cos 0 COS ^)/c \ 

x^l-e J 

0 , 

- sin (P sin C (sin C sin 9cos(p + cos 6> cos C) _ g- 2 ^*/L 2 (i+sinC sine cos ,^+cose cos c)/c 


1 + sin ((sin 9 cos (p + cos 9 cos ( 




(31) 


(32) 

(33) 

(34) 

(35) 

(36) 

(37) 

(38) 

(39) 

(40) 

(41) 

(42) 


The second (exponential) term inside the bracketed term 
at the end of each of these expressions is the contribution 
from the pulsar term. We are in general interested in the 
regime yj = 2'kJLi/c >> 1 (/ = 1,2), and we will present 
results below to leading order in this limit. In the GR 
case, this limit is equivalent to setting the pulsar term 
equal to 0 in the above expressions, i.e., replacing the 
whole bracketed term by 1. This is also the correct thing 
to do for the breathing modes, but more care is needed 
for the other non-GR modes as the term multiplying the 
pulsar term is singular at cos 9 = —1, so we leave this 
term in for now. We will use the above expressions for the 
response functions in Sec. Ill when deriving the overlap 
reduction functions for the different polarisation states. 


III. OVERLAP REDUCTION FUNCTIONS 

The statistical properties of a Gaussian-stationary 
background are encoded in the quadratic expectation 
values of the Fourier components of the waveform, e.g., 
{hA{f,k)h\,{f',k')), where A = {+, x, B, L, X,Y}, in a 
decomposition with respect to the polarisation basis ten¬ 
sors. For an uncorrelated, anisotropic background these 


quadratic expectation values take the form 

{hA{f,k)h%{f,k')) = HA{f)PACk)SAA'S\k,k')6{f-f), 

(43) 

where HA{f) and PA{k) encode the spectral and angu¬ 
lar properties of the gravitational-wave polarisation, 
respectively. [We are assuming here that the spectral 
and angular dependence of the background factorize as 
PA{k)HA{f)] If the background is unpolarised then 
there is the restriction = Px and Px = Py, and 
similarly for iJ+, Hx, and Hx, Hy- 

The functions PA{k) define the anisotropic 
gravitational-wave power distribution on the sky 
for polarisation A, and can be expanded as sums of 
scalar spherical harmonics 

OO I 

PA{k) =Yx Tx Pl^nJlm(k) ■ ( 44 ) 

l—O m— — l 

The expectation value of the correlation between two de¬ 
tectors, labelled 1 and 2, can be written in the form 

/ OO 

( 45 ) 

21 
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where the overlap reduction function, r"^(/), is given by 


r^(/)=E E 


(46) 


with 


l—O m— — l 


^tmU) = 


is^ 




(47) 


Note that a repeated polarisation index A, as in the last 
two equations, is not summed over, unless explicitly indi¬ 
cated with a summation sign. Note also that to simplify 
the notation, we have not included a 12 subscript on the 
overlap reduction functions, as we did in [35], to indicate 
the two pulsars. 

In the following subsections we calculate the overlap 
reduction functions, r^(/), for each mode of the power 
distribution and for each polarisation state, by evaluating 
the right-hand side of Eq. (47) and using the expressions 


for the response functions R^{f,k) given at the end of 
Sec. jn] It turns out that we are able to derive analytic 
expressions for the overlap reduction functions for the -f, 
X, breathing, and two vector-longitudinal polarisation 
modes. For scalar-longitudinal backgrounds, we are able 
to do the (((-integration of (|47|) analytically, but need to 


resort to numerical integration to do the integral over 9. 
Details of the calculations are given in several appendices. 
Plots of as a function of the angle between the 


two pulsars are given in Figs. and We only 

show plots for TO > 0, since = (—l)™r^_^ as a 

consequence of Y/m(fc) = {—l)'^Yi-rn{k)- 


A. Transverse tensor backgrounds 


Analytic expressions for the overlap reduction func¬ 
tions T^if) for uncorrelated, anisotropic (-I-, x) tensor 
backgrounds in GR were derived in |3Sj . For such back¬ 
grounds, we can work in the limit 27r/L/c ^ 1 and set 
the pulsar terms to zero (for which the frequency de¬ 
pendence goes away), obtaining finite expressions for the 
overlap reduction function, even for potentially trouble¬ 
some cases such as cos^ = ±1. Appendix |F| summarizes 
the key analytic expressions derived in that paper. Plots 
of for Z = 0,1, 2, 3 and to > 0 as a function of the an¬ 
gle between the two pulsars are shown in Fig. 0 [rL = o 
as a consequence of (/, fc) = 0 in the computational 
frame.] 


B. Scalar-transverse backgrounds 

For scalar-transverse (breathing mode) backgrounds, 
we can again make the assumption 2'KfL/c ^ 1 and set 
the pulsar term to zero. It then follows that 


— * 
-L Im . 


4 
ttN, 


da; y dcj) (1 — a:) ^1 — a:cosC — \/l — a;^ cos ^ sin A,'”P;'”(a:)e™‘^ 

2(5r„o(l - a;)(l - xcosC)Pi{x) - {5miP^{x) + 6m-iPr^{x))\/l - a;2(l - a;) sin(C 


da: 


= nNrSmO 


1 


1 


1 -f - cos C ) (5/0 - 2 (1 + cos C)<yii + cos C5/2 


-f 7rAj'"'(-l) 2 ' I ^,5,^ - ^5(2 


1 


1 




where we have used the definition of the scalar spherical 
harmonics given in Eq. (All of App. [A] and properties 
of the associated Legendre polynomials summarised in 
App. 0 We see that we are only sensitive to modes of 
the background with I < 2 and |to| < 1. Plots of P;^ for 
Z = 0,1, 2, 3 and to > 0 are shown in Fig. 


not included. We must therefore include the pulsar term 
when evaluating the overlap reduction function for back¬ 
grounds of this form. Using the notation yi = 27r/Li/c, 
?/2 = 27r/L2/c, where Lj is the distance to pulsar /, the 
overlap reduction function for a given (Zto), is given ex¬ 
plicitly by 


C. Scalar-longitudinal backgrounds 


The response for a scalar-longitudinal background, 
Eq. (371, is singular at cos 9 = — 1 if the pulsar term is 

















le angle between the two pulsars ft 
IS are identically zero for Z > 3 or |77 
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(49) 


where 


^ ,, (v/l — sinCcos(/> + xcosC)^ 

Im{v,x) = d(p- 


10 


1 + X cos ( + Vl — x^ sin ^ cos </> 


( 1 _ iy{l-\-x cos — sin cos (p) \ imcf: 

e je 


(50) 


1 

c 

= 0 

C = TT 


Real Imaginary 

Real 

Imaginary 

0 

261 

117 

3.31 

0.254 

1 

-445 

-201 

-6.78 

-0.388 

2 

561 

254 

6.19 

0.567 

3 

-639 

-290 

-6.44 

-0.590 


TABLE I: Values of the co-directional — 0) and anti- 
directional ((■ = tt) overlap reduction function for a scalar- 
longitudinal GW background given for 1 = 0,1, 2, 3. The pul¬ 
sars have y\ = 100 and 2/2 = 200. The values in the table 
correspond to m = 0 modes since all other values of m give 
zero overlap reduction function values. 


The integral for Im{y,x) is challenging to evaluate in 
general; however see App. m for an approximate ex¬ 
pression, valid for large y. As shown in Apps. and 
|lj it can be more simply evaluated for co-directional 
pulsars (i.e., cos^ = 1) and for anti-directional pulsars 
(i.e., cosC = —1). Using the approximate expression for 
Imiy, x) evaluated in App.[^ we then do the integration 
over X given in Eq. (49) numerically. The results of this 
semi-analytic calculation for d'finif) for ^ = 0,1, 2, 3 and 
|m| > 0 are shown in Fig. For these plots we have 
chosen yi = 100 and ?/2 = 200. 


The semi-analytic calculation agrees quite well with 
the full (0, (j)) sky integration, as shown in Fig. (The 
2-dimensional sky integration was actually done using a 
HEALPix |27] pixelisation of the sky.) This plot shows 
the fractional percentage difference between the values of 
the I — 0, m = 0 overlap reduction function roo(/) cal¬ 
culated using these two methods. As can be seen from 
the figure, the agreement is best for values of ( that stay 
away from C = 0 and (^ = n. However, at those spe¬ 
cial points we can use the analytic expressions given in 
Apps. Hand|lj and these are tabulated for / = 0,1, 2,3 in 
Table Q This allows us to obtain a good approximation 
to the overlap reduction function for all (. We note that 
Fig. ID shows that the percentage difference between the 
numerical and semi-analytic curves becomes smaller for 
larger values of yi and y 2 , which is consistent with the 
semi-analytic expression being valid for large y. 


D. Vector-longitudinal backgrounds 


If we ignore the pulsar term, then the response for a 
vector-longitudinal background, Eq. (391, looks singular 
at cos 9 = — 1. However, due to the factor of sin0 in the 
numerator this is a l/Vl -|- cos 9 type singularity which is 
integrable. We can therefore also ignore the pulsar term 
for these backgrounds and obtain a finite result. The 
analytic calculation is very similar to that in App. E of 
[55] for the standard (-I-, x) tensor backgrounds of GR. 
Details of the calculation are given in App.j^ Plots of P^ 
for Z = 0,1, 2,3 and to > 0 are shown in Fig. [P))^ = 0 


as a consequence of RY (/, fc) = 0 in the computational 
frame.] 

We note that in the limit cos^ —>■ 1, the to = 0 
overlap reduction functions diverge. This is because 
in that limit the singularities at (1 -|- fc • ui) = 0 and 
{1 k ■ U 2 ) = 0 coincide and behave like 1/(1 -I- cos0) 
rather than 1 /v/l + cos 9. Again, this singularity is elim¬ 
inated if the pulsar terms are included in the integrand 
and the pulsars are assumed to be at finite distance. De¬ 
tails of that calculation are given in App. |J 1[ 


IV. MAPPING THE BACKGROUND 

In |55| we applied the methodology used to charac¬ 
terise CMB polarisation to describe gravitational-wave 
backgrounds in general relativity. This involved expand¬ 
ing a transverse tensor GR background in terms of (rank- 
2) gradients and curls of spherical harmonics, which are 
closely related to spin-weight ±2 spherical harmonics. As 
described in Sec. |HG| we can use a similar decompo¬ 
sition to represent arbitrary backgrounds with alterna¬ 
tive polarisation states. As explained earlier, for scalar- 
transverse and scalar-longitudinal backgrounds, we ex¬ 
pand in terms of the ordinary (scalar) spherical harmon¬ 
ics, while for vector-longitudinal backgrounds we must 
expand in terms of spin-weight ±1 spherical harmonics. 

In the following subsections, we derive analytic ex¬ 
pressions for the pulsar response functions Rq^)(/) de¬ 
fined in Eq. ( |2^ , for each mode of a background with 
each of the different polarisation states, labeled by P = 
{G,C,B,L,Vg,Vc}- We calculate the response in the 
“cosmic” reference frame, where the angular dependence 
of the gravitational-wave background is to be described. 
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1=0, m=0 



1=0, m=0 




1=1, m=0 
1=1, m=l 




- 1=2, m=0 

radians — — - 1=2, 111=1 

— ■ 1=2, m=2 




- 1=3, m=0 

— 1=3, m=l 

— ■ 1=3, m=2 
—■ 1=3, m=3 


FIG. 3: Plots of the real part (left column) and imaginary part (right column) of for Z = 0,1, 2, 3 as a function of 

the angle between the two pulsars for an uncorrelated, anisotropic background. These were calculated using the semi-analytic 
approximation described in the main text. For these plots we have chosen j/i = 100 and y 2 = 200, where yi = 2nfLi/c and Lj 
is the distance to pulsar I. 


The origin of this frame is at the SSB and a pulsar is 
located in direction it, with angular coordinates (C,x)i 
i.e., 

u°- = (sin cos X, sin C sin x, cos C), (51) 


and is at a distance L from the SSB. In this frame, we 
can again make the approximation x « 0 for the de¬ 
tector locations (i.e., radio receivers on Earth). As was 
done in [35], it is simplest to evaluate the response in 


(51) 
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FIG. 4: Fractional percentage difference between the valnes of the I = 0, m — 0 overlap redaction fnnction roo(/) calculated (i) 
semi-analytically (i.e., using the analytic expression for Iim{y,x) derived in App. and doing the a:-integration numerically), 
and (ii) doing the full (9, (j)) sky integration numerically using a HEALPix |37 | pixelisation of the sky. The dotted curve is for 
yi = 10, t /2 = 20; the dashed curve is for t/i = 50, y 2 ~ 100; and the solid curve is for yi = 100, y 2 = 200, where yi — 2T:fLilc 
and Li is the distance to pulsar 1. Note that the percentage difference decreases as yi and j /2 increase. The vertical dashed grey 
lines at the left and right-hand edges of the plot correspond to the minimum and maximum angular separation (0.95 degrees 
and 174 degrees, respectively) over all pairs of pulsars in the European Pulsar Timing Array (EPTA). 






FIG. 5: Plots of F^ for 1 = 0,1, 2, 3 as a function of the angle between the two pulsars for an uncorrelated, anisotropic 
background. 
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the cosmic frame by making a change of variables of the 
integrand of Eq. ( [^ , so that u points along the z-axis. 
This corresponds to a rotation defined by the Euler an¬ 
gles (a,^, 7 ) = (XjCjO)- Using the transformation prop¬ 
erties of the tensor spherical harmonics y^rn)ab(^) under 
a rotation, it follows that 

Rilm) if) =yim iu)nf {27T fL/c) , (52) 

where TZf‘{2TrfL/c) is proportional to the m = 0 com¬ 
ponent of the response function calculated in the rotated 
frame (with the pulsar directed along the z-axis): 


Note that we need only consider the m = 0 component, 
since the pulsar response must be axi-symmetric in the 
rotated frame, while the tensor spherical harmonics we 
consider are all proportional to in this frame. Thus, 
we see from Eq. (52) that the dependence on the direction 
to the pulsar is given simply by Yimiu), while the distance 
to the pulsar is responsible for the frequency-dependence 
of the response function. Finally, using Eq. (29) with 
X « 0 and doing the integration over </>, we find 


TZ[{2nfL/c) 



(53) 


(27r/L/c) = 2tt 


Att 


dx 


1 1 


2/ -)- 1 J_i 2 \ X 


-Yr 


{ 10)2 


.( 0 , 0 ) ( 


1 — e 


— i27T’fL{l-\-x) j c 


). 


(54) 


where x = cos0. It is this function that we need to 
evaluate in the following subsections. 

We finish this subsection by noting an important result 
implicit in Eq. (52) connected to the distinguishability of 
different background polarisation states. For every po¬ 
larisation type, the response of a pulsar factorises into 
a piece that is dependent on pulsar position, which is 
Yimiu) for all polarisation types, and a piece that de¬ 
pends only on the distance to the pulsar. Even if we had 
infinitely many pulsars distributed across the sky, at any 
given frequency, the best we could do would be to con¬ 
struct a pulsar response map across the sky and decom¬ 
pose it into (scalar) spherical harmonics. The coefficient 
of each term would be a sum of the TZf (27r/L/c)’s for all 
polarisation states, P, which at face value means that it 
would not be possible to disentangle the different polar¬ 
isation states. However, as we will see below, a scalar- 
transverse and transverse tensor background can always 
be distinguished as current PTAs operate in a regime 
in which the response functions are effectively indepen¬ 
dent of the pulsar distance, i.e., the pulsar term can be 
ignored. In that limit, we are only sensitive to modes 
with ^ < 2 of scalar-tensor backgrounds, while trans¬ 
verse tensor backgrounds can only contain modes with 
I > 2. The longitudinal modes cannot be distinguished 
from the transverse modes, however, unless we have sev¬ 
eral pulsars, at different distances, in each direction on 
the sky. For the longitudinal modes the finite-distance 
corrections introduced by the pulsar term are important 


for typical pulsar distances of current PTAs, which gives 
an additional handle to identify those modes. Alterna¬ 
tively, if we made some assumption about how the back¬ 
ground amplitude was correlated at different frequencies, 
e.g., that it followed a power law, we would also break 
this degeneracy as the response of the array to longitudi¬ 
nal modes has a frequency dependence through the same 
term. Thus, it is in principle possible to disentangle every 
component of the background for each polarisation state 
at each frequency, given sufficiently many pulsars at a 
sufficient variety of distances along each line of sight. In 
practice, a pulsar timing array containing Np pulsars can 
only measure 2Np real components of the background at 
any given frequency [33 [SS] and so the resolution of any 
reconstructed map of the background will be limited by 
the size of the pulsar timing array. Roughly speaking, 
to probe an angular scale of the order 1 /Zmax we would 
require Np = (^max -|- 1)^ — 4 pulsars, if we assumed the 
background was consistent with GR and therefore con¬ 
tained only transverse tensor polarisation modes. If we 
allow for arbitrary polarisations we would expect to need 
Np = 3(/niax + 1 )^ pulsars, since we now have structure 
down to / = 0 , and we effectively have three different 
possible polarisation states — transverse (either scalar or 
tensor, but they are distinguished by the I of the mode), 
scalar longitudinal or vector longitudinal. A full inves¬ 
tigation of what can be measured in practice is beyond 
the scope of this current work and we leave it for future 
study. 
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A. Standard transverse tensor backgrounds 


In [35], the standard transverse tensor modes of GR were expanded in terms of gradient and curl tensor spherical 
harmonics, and the corresponding response functions were calculated to be 


« 27r(-l)' « 0, 


(Im) \ 


(55) 


where is a normalisation constant defined in Eq. (C2| of App[^ and the « signs means that the pulsar term 
was ignored for this calculation. Extending the analysis given in |35j to include the pulsar term, we find 


Rtm)U) = YUu)n?i2T^fLlc) , = 0, 


where 


(2) M V / \ 

TZ'riy) = 27T^^ J dx (l-x)(l-x2) 


■ 1 

'-1 


dx^ 


Integrating Eq. (57) by parts twice, 


7^p(y) 


(2 - 2iy + y^)ji{y) - i(6 + iiy + - (6iy - 

ay dy^ dy'^ 


d'^ji . 2 d^ji 


(56) 


(57) 


(58) 


where ji{y) denotes a sp h erical Bessel function, as defined in App. and dj//dy, d'^ji/dy'^, and d^j;/dy^ can be 
simplified using Eqs. (E9)-( ETT] ). Taking the usual limit that the pulsar is many gravitational-wave wavelengths from 
the Earth {y ^ 1), we find TVf\y) ~ 27r(—1)* which is consistent with Eq. (551, where the response functions 
were calculated without the pulsar term. 


B. Scalar-transverse backgrounds 

Repeating the calculation in |35j for an arbitrary scalar-transverse (breathing mode) background, we find 

Rfi^){f) = Yi^{u)nf{2TrfL/c), 


(59) 


with 


nf{y)= 2 ^— 


/ 1 ^(1 - x)Pi{x) (l - e 


= 2^^ 1^,0 - - (-z)'e-^^ 


l-i- ] ji{y) + iji+i{y) 

y 


(60) 


r 


where we used Eqs. (E2), (E9) from App. E to get the 


terms involving the spherical Bessel functions. Since the 
spherical Bessel functions behave like 1/y for large y, the 
terms in square brackets tend to zero as y —)■ oo, leading 
to the approximate expression for the response function 




1 

71 


^10 — 


(61) 


which is valid in the limit where we ignore the pulsar 
term. 

Equation (61) contains a key result of this paper. In 
the limit that y —^ oo, where the influence of the pulsar 
term tends to zero, we find that PTAs will completely 
lack sensitivity to any angular structure beyond ^ = 1 in 
a gravitational-wave background with scalar-transverse 


polarisation. We can verify this analytic result through 
numerical map making and recovery. Using 


- OO / 

hsifM = afm{f)Yim(k). 


(62) 


l—O m— — l 


which relates the expansion coefficients /ib(/,/ c) and 
if) ™ polarisation and spherical harmonic bases 
(see Secs. IIB IIC), we generate a random scalar- 


transverse (breathing mode) background with angular 
structure up to and including I — 10. This injected map 
is shown in the left panel of Fig. To compute the PTA 
response to such a background, we generate a random ar¬ 
ray of Np = 50 pulsars scattered isotropically across the 
sky. We work in the polarisation basis rather than the 
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spherical-harmonic basis here, since the PTA response to 
different angular scales in the GW background is trivial 
in the latter, and we seek a numerical confirmation of 
Eq. (60). The PTA response is computed (with a sky 
resolution set by a given number of pixels iVpix) using 
the Earth term component of Eq. (28), by taking the dot 
product of the array response matrix, R, with the vector 
of amplitude values at each sky-location, h. The matrix 
R has dimensions {Np x A^pix), with each element corre¬ 
sponding to the response of a particular pulsar to gravita¬ 
tional waves propagating in a certain direction (denoted 
by a map pixel), as given by the integrand of Eq. (28). 
The resulting vector is the signal observed by the full 
array, r = Rh. We can invert this in a noiseless map 
recovery by taking the dot product of the Moore-Penrose 
pseudoinverse of R with this observed signal vector. The 
recovered scalar-transverse sky is shown in the right-hand 
panel of Fig. where we note a lack of small-scale angu¬ 
lar structure. We compute the angular power spectrum 
of the recovered and injected maps via HEALPix 
which is capable of rapid map decompositions. The re¬ 
sults are shown in the left-hand panel of Fig. where 
we see that despite the injected map having structure up 
to I = 10, the recovered map only contains structure up 
to and including the dipole. This numerical result is a 
confirmation of the corresponding analytic computation 
inEq. 

We can also check Eqs. (59) and (60), which imply that 
the PTA response to a scalar-transverse background will 
extend beyond the dipole for pulsars at finite distances. 
We do so again with numerical map making and recovery, 
by using the full Earth and pulsar term scalar-transverse 
response function given in Eq. ( |28[ ). The pulsar term will 
be highly oscillatory across the sky, so we expect some 
numerical fluctuations in our results. For this study we 
inject white Gaussian noise in each pulsar measurement, 
with an amplitude such that the GW background remains 
in the strong signal limit. In the right-hand panel of 
Fig.0 we see that the PTA has increasing sensitivity to 


higher multiple moments in the background as y is in¬ 
creased. At y ~ 5 — 10 the PTA is able to recover the 
full angular structure of the background, but also suffers 
from noise leakage at higher multipoles, since the non¬ 
zero response of the pulsar term at these higher multi¬ 
poles amplifies noise arising from the pixelation of the 
sky. The pulsar term response peaks at Z ^ y, such that 
for PTAs with y = 15, 20 we see a drop-off in sensitivity 
at I ~ 15, 20, even though the response is merely amplify¬ 
ing pixel noise at these multipoles. For y ^ 20 the Earth 
term behaviour is recovered, and we observe a lack of 
sensitivity to modes beyond dipole. To put these results 
into context, we recall that y = 27r/L/c and peak PTA 
sensitivity to a gravitational-wave background occurs at 
/ ~ 1/T where T is the total observation time. For 
T = 20 years, this gives / ^ 1.6 nHz. Thus in order for a 
PTA to have sensitivity to structure in a scalar-transverse 
sky beyond dipole, we need y < 10, which corresponds to 
all pulsars in our array being at a distance of < 0.01 kpc 
from Earth. Given that most timed millisecond pulsars 
have distances > 0.2 kpc, it is unlikely that this extended 
reach to sensitivity beyond dipole modes will be possible 
with current arrays. 

Using the mapping response functions Rp,„)(/) calcu¬ 
lated above, we can also compute the overlap reduction 
function for an uncorrelated, anisotropic background, re¬ 
covering the result given in Sec. mE Details of that 
calculation are given in App. 


C. Scalar-longitudinal backgrounds 

For an arbitrary scalar-longitudinal background we 
find 


Rfir^){f)=YUu)nfi27rfL/c), (63) 


where 


where Hi{y) 
terms in the 


^f(j/) = 27r (l-e 

^lix) (l - 


'-1 

/■' 1 
= 27r / da; - 


-1 -I- a: -I- 


1 -I- a; 


(64) 


= 27^<j-5^o+3<5^l + (-^)'e-*^ 


1 - *- ) jiiy) + iji+i{y) 


Imy) 


is defined by Eq. (H6) in App. Since the spherical Bessel functions behave like 1/y for large y, the 
square brackets above tend to zero as y —>■ oo, yielding 




-Sio + -f ^Hi{y) 


(65) 


This is valid for y ^ 1, but y finite. 
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FIG. 6: Maps of the real amplitude component of a scalar transverse (breathing mode) background. (Left) a randomly generated 
scalar-transverse gravitational-wave sky, with structure up to and including I = 10. (Right) the corresponding recovered sky, 
computed by first forming the observed signal vector for an array of Np = 50 pulsars via r = Rh, where each element of the 
array response matrix, R, corresponds to the response of a particular pulsar to gravitational waves propagating in a given sky 
direction. We perform a noiseless map recovery by computing R+r (where R'*' is the pseudo-inverse of R) which gives the map 
in the right panel. We note the lack of small-scale angular structure in the recovered map compared to the injected map. 




(a) 


(b) 


FIG. 7: (Left) A comparison of the angular power spectra of the injected scalar-transverse sky map shown in the left-hand panel 
of Fig.[^ and the PTA-recovered map shown in the right-hand panel of the same figure. We see that PTAs will completely lack 
sensitivity to angular structure in a scalar-transverse gravitational-wave sky beyond the dipole level. This result is confirmed 
analytically in Eq. (|61||. (Right) We use the full Earth and pulsar term response from Eq. (281 to investigate map recovery with 


finite y. The pulsar-term will be highly oscillatory across the sky, so we expect some numerical fluctuations in our results. As y 
is increased the PTA shows greater sensitivity to higher multipole moments in the GW background. At t/ = 10 the PTA is able 
to recover all modes in the injected map, although the non-zero sensitivity of the pulsar-term response at higher multipoles 
amplifies noise from the pixelation of the sky. For y > 20 the Earth term behavior is recovered, and we observe a lack of 
sensitivity to modes beyond dipole. See text for further details. 


D. Vector-longitudinal backgrounds 


As discussed in Sec. |II B[ we can expand each Fourier component of a vector-longitudinal background in terms of 
gradient and curl tensor spherical harmonics ’''^^lich are simply related to the spin-weight ±1 
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spherical harmonics defined in App. It is convenient to relate this expansion 


OO L 

ha,{f,k)=Y. E 


1—1 m— — l 

to a similar expansion in terms of the polarisation basis: 

KbU, k) = hxif, k)e^i,{k) + hrif, k)e^i,ik) ■ 

The relationship is 

1 

I 

Im 


hxif, k) ± zhrif, k) = T^J2 (<-)(/) ± *« m (/)) ±iYimik), 


( 66 ) 


(67) 


( 68 ) 


or, equivalently, 


hxif,k) = 2 ^E Kh(/) - i^Uk) 


-.Vc 

{Im) 


if) (^-lYlmik) + lYimik) 

hvif, fc) = A E {-lYiM - lYimik)) + [-^YUk) + lY^k) 


(69) 


where iiTzm)^) are the spin-weight ±1 spherical harmonics dehned in App. 

The expressions for the grad and curl response functions for an arbitrary vector-longitudinal background can be 
calculated using the same methods as in the preceding subsections. We find 


where 


Rnim)U) = YUuiin^^i^xfLjic ), <((,„)(/) = o, 


TZY‘^ iy) = J a;(l — a:) ^1 — e 


(70) 


(71) 


Thus, the response to vector curl modes is identically zero for pulsar timing arrays, as is the case for tensor curl 
modes, as shown in [^. Evaluating the integral in Eq. 0 by parts we find 




-2(5;o + ^<5(1 + (-l)'e Ax il + i2 + iy)x + Piix) 

il + l)jiiy) -iy- ii2l + 3))ji+i(?/) - iyji+ 2 iy) 


= 7T^^^Ni\-dn+2i-iye-^y 


1 - 




(72) 

(73) 


where we have dropped the Sig term since for spin-weight 
±1 harmonics we have ^ > 1. Taking the usual limit that 
the pulsar is many gravitational-wave wavelengths from 
the Earth, y ^ 1, and using the asymptotic result 


ji(y) » Ein -1-0 , fory»l, (74) 

we find 


<^(/) « 27rYUu) 


^6n + (-1)' ('H 


(75) 


E. Overlap reduction function for statistically 
isotropic backgrounds 


For a statistically isotropic, unpolarized and parity- 
invariant background (see, for example, Eqs. (52)-(54) 
of [55] 1 


Tif) =J2CiTiif), (76) 

i 


where 


As expected, this agrees with the result obtained by eval¬ 
uating the integral in Eq. (71) without the pulsar term, 
i.e., making the replacemeiitj 1 — exp[—iy(l -|- x)]} —>■ 1. 


r*(/)= E J2^Hlm)if)R2ilm)if)- ( 77 ) 

m— — l P 
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Here is a sum over the polarization states for a 
particular type of background (e.g., P = {Vq^Vc} or 
P = {G,C} for vector-longitudinal or transverse ten¬ 
sor backgrounds). Using the results of the previous sub¬ 
sections, we have in the limit j/i S> 1, j /2 ^ 1 (where 
yj = 2TrfLj/c as before): 

Transverse tensor modes {I > 2): 


which was found in |35j . 
Scalar-transverse mode {I > 0): 


rf(/)«7r(2; + i)i 


^10 + glJ/l 


Pi (cos C) ■ 


rf (/) « 7r{2l + l){Nf)'^Pi{cos Q , (78) Scalar-longitudinal mode {I > 0): 


(79) 




7r(2Z-h 1 ) 



1 - ^ (Hoiyi) + H^{y2)) 


+ Sii 


Vector-longitudinal modes {I > 1): 






Pi (cos C). (81) 


Note that only the scalar-longitudinal overlap reduction 
functions Ff (/) are actually frequency-dependent in the 
large y limit, via their dependence on P[i{yj). The other 
overlap reduction functions depend only on the angular 
separation C, between the pair of pulsars. 

As shown in [35], an isotropic, unpolarized and uncor¬ 
related background has Ci = 1 for all 1. In Fig. we plot 
approximations to T®, T^, T^, and T^ correspon din g 
to different values of Z^ax in the summation of Eq.(76l, 
taking Ci = 1 for all I up to Imax- (Recall that for the 
vector overlap reduction function, the summation starts 
at Z = 1, while for the tensor overlap reduction function, 
it starts at Z = 2.) These hnite Zmax expressions are com¬ 
pared to the Z = 0, TO = 0 components of the overlap 


reduction functions calculated in Sec. Ill and plotted in 
Figs. BHii The normalization is different than in 
those hgures, since the Z = 0, to = 0 components need to 
be multiplied by •\/47r/2 in order to obtain the isotropic 
overlap reduction function. (The factor of \/47r comes 
from Foo(^) = l/VdTr; the factor of 1/2 is needed to 
get agreement between Eq. (43) and Eq. (32) of |3S| for 
isotropic, unpolarized backgrounds.) 

Figure confirms what was found for the transverse 
tensor modes in |35j . namely that a good approximation 
to the full overlap reduction function can be obtained by 
including only a relatively small number of modes in the 
sum. The maximum Z required in the sum is approx¬ 
imately 1,4,10 and 20 for the scalar-transverse, trans¬ 
verse tensor, vector-longitudinal and scalar-longitudinal 
backgrounds respectively. 


V. CONCLUSION 

In this paper we have investigated the overlap re¬ 
duction functions and response functions of PTAs for 


l + lmyi) + ff^(y 2 )) 


r 


+ ^^!(yi)ffi*(y2)j^i(cos(). 

(M 


non-GR polarisations of gravitational waves. The over¬ 
lap reduction function describes the sensitivity of a 
pair of pulsars to a gravitational-wave background in a 
cross-correlation analysis. The cross-correlation signa¬ 
ture traced out by the overlap reduction function from 
an entire array of precisely-timed millisecond pulsars 
will aid in isolating any gravitational-wave signal from 
other stochastic processes which may have similar spec¬ 
tral properties. Hence, current searches for stochastic 
gravitational-wave backgrounds rely on models of the 
overlap reduction function as the smoking-gun signa¬ 
ture of a signal. For an isotropic stochastic background 
in GR, the overlap reduction function is known as the 
Hellings and Downs curve, and depends only on the angu¬ 
lar separation between pulsars in the array. The overlap 
reduction functions for arbitrary anisotropic stochastic 
backgrounds in GR were investigated in Mingarelli et al. 
[30] , Gair et al. [35] , where it was shown that these func¬ 
tions are now dependent on the positions of each pulsar 
relative to the distribution of gravitational-wave power 
on the sky. 

The gravitational wave polarisation has a strong in¬ 
fluence on the overlap reduction function through the 
form of the pulsar response functions. Ghamberlin and 
Siemens |34| studied the form of the overlap reduction 
functions for isotropic backgrounds of gravitational waves 
for scalar-transverse, scalar-longitudinal, and vector- 
longitudinal polarisation modes. In this paper, we have 
extended that analysis to find analytic expressions for the 
overlap reduction functions for anisotropic non-GR back¬ 
grounds. A key result of this work is that PTAs will com¬ 
pletely lack sensitivity to structure beyond quadrupole 
in the power of a scalar-transverse background. This re¬ 
sult holds regardless of the number of pulsars, timing- 
precision, or observational schedules—it is a property of 
the geometric sensitivity of PTAs to gravitational-wave 
signals of scalar-transverse polarisation. Additionally, we 
have found analytic expressions for the overlap reduction 
functions for arbitrary anisotropic vector-longitudinal 
backgrounds. We also derived a semi-analytic expression 
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— Imax = 1 

— Imax = 2 

— Imax = 3 

— Imax = 4 

— Imax = 5 

— Imax = 10 

— Imax = 20 

— 1=0, m=0 


FIG. 8: Approximations to the overlap reduction functions for an isotropic, unpolarized and uncorrelated stochastic background, 
plotted as a function of the angle between a pair of pulsars. The approximations are obtained by summing products of the 
response functions over I for different values of Zmax. Panel (a): transverse tensor background; Panel (b): scalar-transverse 
(breathing) background; panel (c): scalar-longitudinal background; panel (d): vector-longitudinal background. We are working 
in the large y limit for all of these cases. For the scalar-longitudinal background, we have taken j/i = 100 and y 2 = 200. The 
thick black line in each plot is the “full” expression for the overlap reduction function, corresponding to the limit Zmax —^ cio. 
(These limiting expressions equal \/47r/2 times the I — 0, m — 0 component of the overlap reduction functions calculated in 
Sec. Ill ) For the scalar-longitudinal case, the full expression was calculated numerically. 


for the overlap reduction functions of anisotropic scalar- 
longitudinal backgrounds, in which case a consideration 
of the pulsar-term is crucial to avoid divergences. 

In the second half of this paper, we extended the for¬ 
malism of our previous work in Gair et al. |35| . where the 
Fourier amplitudes in a plane-wave expansion of the GR 
metric perturbation were decomposed with respect to a 
basis of gradient and curl spherical harmonics, which are 
related to spin-weight ±2 spherical harmonics. By deter¬ 
mining the components of the background in such a de¬ 
composition it is possible to construct a map of both the 
amplitude and the phase of the gravitational wave back¬ 
ground across the sky, rather than simply reconstructing 
the power distribution. The decomposition in terms of 
spin-weight ±2 spherical harmonics is made possible by 
the transverse-traceless nature of the GR gravitational- 
wave metric perturbations. Here we have appealed to 
the structure of the gravitational-wave metric pertur¬ 
bations for non-GR polarisations to perform the same 
procedure—the Fourier amplitudes of scalar modes can 
be expanded in terms of ordinary spin-weight 0 spherical 
harmonics, while the vector mode amplitudes can be ex¬ 


panded in terms of a spin-weight ±1 spherical harmonic 
basis. In so doing, we found that PTAs lack sensitivity 
to structure in the polarisation amplitude of a scalar- 
transverse background beyond dipole anisotropy, which 
can be used to explain the lack of sensitivity to power 
anisotropies beyond quadrupole. This result was veri¬ 
fied through numerical map making and recovery, where 
we found some sensitivity to modes beyond dipole when 
y = 2'KfL/c was very small, but this would require all 
pulsars to lie within a distance of 0.01 kpc from Earth. 
We also found that PTAs will lack sensitivity to vector 
curl modes for a vector-longitudinal background, which is 
analogous to the finding in Gair et al. [35] that PTAs are 
insensitive to the tensor curl modes of gravitational-wave 
backgrounds in GR. 

This paper provides several ready-to-use expressions 
for overlap reduction functions for non-GR stochastic 
backgrounds with arbitrary anisotropy. These expres¬ 
sions can be trivially plugged into any current or planned 
PTA stochastic background search pipeline to obtain lim¬ 
its on the strain amplitude of a non-GR gravitational- 
wave sky. We also provide several ready-to-use expres- 





















19 


sions for the response functions of a single pulsar to 
anisotropies in a non-GR gravitational-wave background. 
The implications of this are that we can use an array of 
pulsars to perform a Bayesian or frequentist search for 
the angular dependence of the Fourier modes of a plane- 
wave expansion of the gravitational-wave metric pertur¬ 
bations, and in so doing produce maps of the polarisation 
content of the sky that include phase information rather 
than simply map the distribution of power. 

The results in this paper also indicate what is possi¬ 
ble to measure in principle with a sufficiently extensive 
pulsar timing array. The dependence of the response on 
the pulsar location on the sky is proportional to Yim(u), 
where ii is the direction to the pulsar, for all polarisa¬ 
tion types. By decomposing the pulsar response map, 
at a particular frequency, into regular (scalar) spheri¬ 
cal harmonics, the coefficients of each Yim{u) mode of 
the response map can be determined, but these coeffi¬ 
cients will be a sum of the contributions from each of 
the polarisation types. Scalar-transverse and transverse 
tensor backgrounds can be distinguished because PTAs 
typically operate in a regime in which the pulsar term 
is negligible and so the response is independent of the 
distance to the pulsar. In that regime, PTAs are only 
sensitive to modes of the scalar-transverse background 
with I < 2, while transverse tensor backgrounds can only 
contain modes with I > 2. However, longitudinal back¬ 
grounds can only be distinguished from transverse back¬ 
grounds if there are multiple pulsars along a given line 
of sight, or if there is a known correlation (e.g., a power 
law) between the background amplitudes at different fre¬ 
quencies. In either of these scenarios, we can exploit 
the dependence of the pulsar term on 27r/L/c, which is 
much more significant for the longitudinal modes of the 
background. Thus, in the limit of infinitely many pul¬ 
sars distributed across the sky at a range of distances, 
we would expect to be able to measure the entire con¬ 
tent of the background in each polarisation state and at 
each frequency. In practice, of course, a pulsar timing 
array of Np pulsars can only measure 2Np real compo¬ 


nents of the background [33 |3B] , and so the resolution of 
any map that we produce will be limited by the number 
of pulsars in the array. Roughly speaking, to produce 
a map of the gravitational wave sky in all polarisation 
states to an angular resolution of 1/Giax would require 
Ap = 3{lmax + 1)^ pulsars, but this should be explored 
more carefully in the future. 

For a further discussion of the prospects of this type of 
mapping analysis, see Gair et al. |35j . Cornish and van 
Haasteren [35] . We plan to apply the results of this paper 
to the analysis of real data, to map the amplitude and 
phase content of non-GR gravitational-wave backgrounds 
influencing the arrival times of millisecond pulsars. This 
will allow us to place constraints on beyond-GR polari¬ 
sations of nanohertz gravitational waves. 
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Appendix A: Spin-weighted spherical harmonics 


This appendix summarizes some useful relations involving spin-weighted and ordinary spherical harmonics, sYim{k) 
and Yim(k). For more details, see e.g., Goldberg et al. [39] and del Castillo [40]. Note that we use a slightly different 
normalization convention than in Goldberg et al. [33|- Namely, we put the Condon-Shortley factor (—I)™ in the 
definition of the associated Legendre functions Pl^{x), and thus do not explicitly include it in the definition of the 
spherical harmonics. Also, for our analysis, we can restrict attention to spin-weighted spherical harmonics having 
integral spin weight s, even though spin-weighted spherical harmonics with half-integral spin weight do exist. 

Ordinary spherical harmonics: 


Yi^{k) = Y^e, <t)) = Nr prices , where 


(2/ -I- 1 (Z — to)! 

47r (Z -|- to)! 


(Al) 
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Relation of spin-weighted spherical harmonics to ordinary spherical harmonics: 


sYlm{e,cj,) = 


sYi^{e,<p) = 


{l + s)l 

l{i + sy. 
{i-sy 


d’^Yim{ 0 , (j)) for 0<s<l, 


(-iyd-^Yira{e,^) for -l<s<0, 


where 


■ O ~ 

dri = -{sm9y -^+icsc6— (sin 6 »)“®? 7 , 
\_ou acp] 

_ ' f) C\ ' 

9?7 = —(sin 0 )“® ——icscd— (sin 0 )®? 7 , 
yaO d(p\ 


and rj = ri{ 0 , 0) is a spin-s scalar held. 


Series representation: 




(I + my.{l — my. 21 - 1-1 . 

{I + 3)1(1 — s)! Itt 


n ^ ^ 


I — s\ / I + 3 
k I \k + 3 — m 


(-l)'-'=-"e*™'^(cot6»/2) 


2fc+s—m 


Complex conjugate: 


Relation to Wigner rotation matrices: 




D^m'ra(y},0,y) = (-1)^ 


21 + 1 




[D‘^,^{ye,y)]* = (-ir 


21 -h 1 


-mYl,m'(0A)e^^^. 


Parity transformation: 


Orthonormality (for hxed s): 


sYirail^ -B,(j) + TT) = (-1)' (6», 0). 


(‘‘2.7^ nTT 


sYim(k) sYy.^,(k) = / d0 / sin6*d6» 0) 0) = ( 5 h/( 5 „ 

Jo Jo 


Addition theorem for spin-weighted spherical harmonics: 


^ sYUe^,<j,y 02 ) = s'YUOs, 


where 


cos 6*3 = cos 6*1 cos 02 + sin 9 i sin 02 cos(02 — 0i), 


^_*(03+X3)/2 _ cos ^(02 - 0i)cos|(02 - 0i) - zsin|(02 - 0i)cos^(0i + 92 ) 

^Jcos^ i (02 - 0 l) COS^ i (02 - 0 l) -I- sin^ ^(02 - 0 i) cos^ i( 0 i -I- 02 ) 
g*03-X3)/2 _ cos ^ (02 - 01) sin ^ (02 - 01) + i sin |(02 - 0i) sin | (0i -H 02) 

-i/cos^ i(02 - 0 i) sin^ i(02 - 0 i) -I- sin^ ^(02 - 0 i) sin^ ^(01 -I- 02) 



Addition theorem for ordinary spherical harmonics: 


21 


‘ 2 / + 1 
^ YUki)YL{h) = P/(fci • fe). 


1 — — 1 


47r 


Integral of a product of spin-weighted spherical harmonics: 

^{2li + 1)(2Z2 + 1)(2/3 -I- 1) f li I 2 I 3 




where 


d (fc) 82^2^3 (^) S3^3”l3 (^) 


47r 


TOl m 2 m 3 


h h h 

mi m 2 m 3 


is a Wigner 3-j symbol, which can be written as 


(A13) 


h h h \ 

—Si —S 2 —S 3 j ’ 

(A14) 


f I I' L \ Ul + l' - L)\{1 -l' + L)\{-1 + 1' + L)\{1 + my.il - m)!(l' + + M)!(L - M)! ^ 

\m m' M J ~ \l + y + ^ 

X f A 1 

zlil + 1' — L — z)!(Z — TO — z)!(Z' + to' — z)\iL — l'-\-m + z)!(L — Z — to' + z)! ’ 

See, for example, Wigner [IT], Messiah [35], Landau and Lifshitz [53] and references therein. Note that although this 
sum is over all integers it contains only a finite number of non-zero terms since the factorial of a negative number is 
defined to be infinite. 


Appendix B: Gradient and curl rank-1 (vector) spherical harmonics 

The gradient and curl rank-1 (vector) spherical harmonics are dehned for Z > 1 by 

r,a.,. = = lew, ("It^^■ 


de 


where 9 and (j) are the standard unit vectors tangent to the 2-sphere 

9 = cos 9 cos (j)x + cos 0 sin (/> y — sin 0 z 
4 > = — sin (px + cos (j) y , 

is a normalisation constant 


Y (Z + 1)! ’ 


and €ab is the Levi-Civita anti-symmetric tensor 


Eab = Q ^ , 9 = detigab) ■ 


(Bl) 


(B2) 


(B3) 


(B4) 


Following standard practice, we use the metric tensor on the 2-sphere gab and its inverse to “lower” and 
tensor indices—e.g., e^b = g'^’^^ab- In standard spherical coordinates ( 0 , (/i), 

(S sin%) ’ = 

The grad and curl spherical harmonics are related to the spin-weight zbl spherical harmonics 




'(Z-l)! Nr b.dPl 


,_ ,±il-xr—^+mPrix)]e^^‘^, wherea; = cos0 

{l + iyVl^Yr^ \ ' dx ' ^ ^ ' 


(B5) 


(B 6 ) 
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via 


y^m)a{k) ± = ±^( 0 a ± ^k) ^lYlr. 


ik) 


or, equivalently, 


y{lm)a{k) 2^/2 
y{lm)aik) — 


(^-lYi^ik) - iri^(fc)) Oa + I [-lYir^ik) + lYim 
(^-lYim{k) - lYiraik)^ k “ * (^-lYi^{k) + lYim 


For decompositions of vector-longitudinal backgrounds, as discussed in the main text, it will be 
construct rank-2 tensor fields 


— y(?m)akb + Yki)bka , 


wFs _ 

^ (lm)ab ^ (lm)a 


~ y(lm)ak + y(lm)bk , 


_ 

^ {lm)ab ^ (lm)a 


where k is the unit radial vector orthogonal to the surface of the 2-sphere: 

k = sin 6 cos (j)X + sin 6 sin ipy + cos 6 z . 
These fields satisfy the following orthonormality relations 

d^f]^ Y^-^^kk)Y^kk’’*Ck) = SwSmm' 


^^^kY^kkk)Ykrnk"k)= 0 - 


Vc ab * / 


is^ 


Appendix C: Gradient and curl rank-2 (tensor) spherical harmonics 

The gradient and curl rank-2 (tensor) spherical harmonics are defined for I > 2 by: 

1 


Y(^lm)ab kYl I n9abY(^im);c 


N, 


Y[lm)ab 2 {Y(^lm)-,ac£ b + Yfj^y^bc^ a) > 

where a semicolon denotes covariant derivative on the 2-sphere, and is a normalisation constant 


y (! + 2 )! ■ 


Using the standard polarisation tensors on the 2-sphere: 

kbik) = daOb - kk : 
kb(k) = kk + kk : 


where 0, (j) are given by Eq. (B2), we have 


„ , (2)yV; r ^ 

Y(fm)abik) = [W(i^)ik)e+^{k) + X(^i^){k)e^^{k) 

„ , (2)yV; r ^ . 

Yk)abik) = W^,^)ik)e:,ik) - Ap^)(A:)e+(fc) 


(B7) 


(B8) 

convenient to 

(B9) 


(BIO) 


(Bll) 


(Cl) 


(C2) 


(C3) 


(C4) 
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where 


W(^lm){k) = ~ ■'■ gjjj2 ^ j Ylm{k) = + l)j Ylm{k) ; 


^ OlTD / Ci \ ^ 

These functions enter the expression for the spin-weight ±2 spherical harmonics 

(2)Ar, r . ^ . 1 

±2Ylm{k) = ^ W(^im){k) ±iX(^im)ik) , 

which are related to the grad and curl spherical harmonics via 

i"(L)a6(fc) ± *>^(L)a6(fc) = ^ (dW ± <,ik)) ^2YUk) ■ 

Note that the grad and curl spherical harmonics satisfy the orthonormality relations 


Js^ 


JS^ 


I 


Y^^y,{k)Y^.^,y^*{k) = 5u'5rnm' 


Appendix D: Legendre polynomials and associated Legendre functions 

The following is a list of some useful identities involving Legendre polynomials Pi{x) and 
functions P™(a;). For additional properties see e.g., Abramowitz and Stegun [15] . 


Differential equation: 




da; 




(l-x 2 ) 


pr{x) = Q. 


Useful recurrence relations: 


^ + ^)pr-i{x) i)pr+,[x)\ , 

y/l-x^^P^{x) = i [(Z + m)(l - m + l)Pr”^(^) - Pr^Hx)] ■ 


Orthogonality relation (for fixed m): 


da;Pr(a;)P,U(a;) = 


J dx Pi{x)Pi>{x) = 


{2l + l){l-m)l 

^ s.,. 


{21 + 1 ) 

Relation to ordinary Legendre polynomials, for m = 0,1, • • • , ^ 

jm 

pr(a;) = (-ir(i - p^x), 


Rodrigues’ formula for Pi{x): 


p-(x) = (-ir|^pr(x). 




(C5) 

(C6) 

(C7) 

(C8) 

associated Legendre 

(Dl) 

(D2) 

(D3) 

(D4) 

(D5) 
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Series representation of Legendre polynomials: 




fc =0 


{I + k)l f 1 + X 


Useful recurrence relation: 


which iterated yields 


(fc!)2(Z-fc)!V 2 J ^ ik\)Hl-k)\ 

{21 + l)xPi{x) = (/ + l)Pi+i(a:) + lPi_i{x), 


(D6) 

(D7) 

(D8) 


^ {l + 2){l + l) „ 4;3 + 6P-1 

^ {21 + 3){2l + ^ {21 + 3)i2l + 1){21 - ^ 4P _ 1 ' 

Appendix E: Bessel functions 


The following is a list of some useful identities involving Bessel functions and spherical Bessel functions of the first 
kind, and ji{y). For additional properties, see e.g., Abramowitz and Stegun |46) . 


Integral representation of ordinary Bessel functions: 


i{n4i+y cos cj>) 


Integral representation of spherical Bessel functions: 


Relationship between ordinary and spherical Bessel functions: 


Plane wave expansion: 


functions: 

2 (-*)'i;(y) = y dx Pi{x)e~’-y^ . 

3ssel functions: 

oo 

-^ 2 ^fk■S/c ^ ^-^ycose ^ + l)P,(cos6>). 


i =0 


Asymptotic behaviour: 


Jn{y) " 

1 


yy 

' r(n+ 1 ) V 

2/ ’ 

Jn{y) " 

V Try 

cos 

/ niT 

[y-^- 

ji{y) - 

1 . 

s - sm 

y 

b- 



O 


for 2 / > 1 , 


A useful recurrence relation: 


Another useful recurrence relation: 


which iterated once yields 


ji-i{y) + ji+i{f) = ■ 

y 


= -jiiy) - ji+i{y ), 

dy y 


d^ji ~ t \ 21 + 1 . . 

= -^^Jiiy) - —Ji+i{y) + Ji+2{y), 


and twice yields 


d^ji l{l-l){l-2) ^ 3P . ^^^3{l + l). 

— =-- ■Ji{y) - -:;;^Ji+i[y) + —y— ji+ 2 {y) - ji+3{y) ■ 


dy^ 


yO 


y" 


y 


(El) 

(E2) 

(E3) 

(E4) 

(E5) 

(E6) 

(E7) 

(E8) 

(E9) 

(ElO) 

(Ell) 



















25 


Appendix F: Analytic calculation of the overlap reduction functions for transverse tensor backgrounds 

For completeness, we include here expressions for the overlap reduction functions for anisotropic, uncorrelated 
backgrounds having the standard transverse tensor polarization modes of GR. These were derived in App. E of |35) . 
Here we present only the final results; readers should consult [33] for details. 

For all I, m: 


rL(/) = o, 


which trivially follows from the fact that (/, fc) = 0 in the computational frame. 


(FI) 


For m = 0: 


^w(f ) = + l)7r| (^1 + ^ cosC^ Sio - i (1 + cosC) Sn + ^cosCfe 

- (1 + cosC).F-o ; o(cosC) - (1 - cosC).F+^,; o(cosC)| • (F2) 

For TO = 1: 


r+(/) = i^/(2ZTlV^/|^|2sinC 


^ i <5,2 

3 5 


(1 + cosC)^^^ T-_ 


(1 — cosC)^^^ 




For TO = 2,3, • • •: 

1 


^ 1 /j— —/ (Z — to)! I (1 + cosC)" , (l + cosO™ T , 

r,„.(/) = -A<^' + i)'V((T^ ^ 


(1 -cosQ ^ .^+ 


^•^m+l,l./,m(cOsC)- 


(1-cosC)™ .^+ 


For TO < 0: 


(1 + cosC) 2 


rL(/) = (-i)'"r+ (/)■ 


(1 + COS^) 2 


2^ —I 771,0,1,771 


The functions which appear in the above equations are defined by 


T. 


q,r,1,771 


(cosC) = 


Kr,LmicOSC) = 


/ dx — -—-— Piix). 

/_i {i-xyax-^ ^ 

r , (i-a;)« d™ , 

/ ax 7"-^-;--tMX) . 

/-cosC y + xYdx-- 


(F3) 


(cosC)>. (F4) 


(F5) 


(F6) 


These functions also arise when calculating the overlap reduction functions for the vector-longitudinal polarization 
modes. The integrals can be evaluated analytically as shown in App. [^of this paper (or in App. E of [35]). 


Appendix G: Evaluating the Im{y,x) 


integral for the overlap reduction function for scalar-longitudinal 
bacgkrounds 


The response for a scalar-longitudinal background, Eq. (37), is singular at cosd = — 1 if the pulsar term is not 


included. We must therefore include the pulsar term when evaluating the overlap reduction function for backgrounds 
of this form. We use the notation yi = 27r/Li/c, 1/2 = 27r/L2/<3, where Lj is the distance to pulsar /, that was 
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introduced in the main body of this paper. I n the following, we will ensure that we keep all terms up to constant 
order (yi)°, {y 2 )^■ The final expression, Eq. (G12), contains some terms of higher order, but these are incomplete. 


This will be discussed further below. The components of the overlap reduction function are given by 


rL(/) = 


dx 


1 — e 


-iyi(l+a:) 


Im{y2,x) 


pr(x).. 


(Gl) 


where 

Im{y x) = [ d(f> ~ ^ sinC cos ^ + a: cos C) (l — e»y(l+2: cos Q+yi-x^ sin C cos 0)'\ ^im<p 

Jo 1+ X cos C + \/l — x"^ sin C, cos 4> ^ ’ 

The integral for Im{y,x) can be simplified by writing 


(G2) 


Im{y,x) = d(^ 


_ gjy(l+a: cos C+\/l-a;2 sin C cos 0)^ 


X COS C + \/l — sin ( cos (j) — 1 
2tt (xcosC - 1) (dmO - Jm{y sin C\/l - 

+ 7rsinC\/l - x"^ Jm+i(ysinC\/l - a;^) - J„_i(ysinC\/l - 


Zl/(l+XCOS C) 


+ im{y,x) 


where 


im{y,x) = d(j) 


2 ^ j ^ Qiy{l-\-x cos (^+-\/l —a;2 sin cos (p) \ 

r\^ ^ _Z, 


0 1+ X cos C + vT^-a? sin C, cos (j) 


im(j) 


(G3) 


(G4) 


and Jn{y) denotes the Bessel function of the first kind. For large values of y, Bessel functions have the asymptotic 
form given in Eq. (E6), and we will use this to drop certain terms when we take the limit y/ —>■ oo later. 

To evaluate the integral Im{y,x), we first note that lm(0,x) = 0 and 


dl 


p27r 

_ I gZm(/i+iy(l+a:cos — sin C cos </)) 

dy Jo 

= J^( 2 /sinC\/l - x^). 


(G5) 


This last equation can be integrated as follows. For 1 + xcos^ ^ sin^Vl — x'^ (which corresponds to x + cos(( ^ 0) 
the integral to infinity can be computed as 


ijn{oo,x) = 27r(-l)’' 


sin 


CVl — 


|m| 


COS C + a;| I 1 + X cosC + |x + cos C| 


(G6) 


This is divergent at x = — cosC, but that is an artefact of taking the limit y —>■ oo. To evaluate Im{y,x) for finite y 
we can write 


pOO 

ini{y,x) = i^{oo,x) + 27:1^"^^ / dy sinCVl - x^). 

J V 


(G7) 


For the range of y in the integral, we can approximate the Bessel function using Eq. (E6). The corrections to this 
approximation take the form of trigonometric functions times factors of and will contribute terms of order 

f/^y and smaller to the result. To obtain a result accurate to at least 0(y?,y2), we therefore just need to evaluate 


TT sm 


CVl — x^ 


dy 

Vy 


iy{l-\-x cos (^) 


COS 


(^y sin C yr 


mn TT 


in CVl — - 


TT sin 


^mgZ7r/4 

+ X — 


) 


Pc 2/(1 + 2;-)^ + 


i/l + x+ 


[Vyi^ + x+) 


(G8) 
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where x± is shorthand notation for 


x± = X cos C ± sin C \/\ — x^ , 


and 


r 

M = / 

'Jy 


FAy)= I due*“'= 


C\ \ -y\+iS 


Here C{x) and S{x) are the Fresnel cosine and sine integrals, defined by 


C{y)= J ducos(^^u^) , S{y) = J dwsin(^^u^^ . 


Thus, 


(G9) 


(GIO) 


(Gil) 


im{y,x) = 27r(-l)’' 


1 


sin 


CVi — 


\m\ 


I COS C + a;| \ 1 + a; cos ^ + |a: + cos C| 


in CVl — : 


TTSin 


\/(l + X-) 


Fr 




\/(l + a^+) 


(v/j/(l + a:+) 


• (G12) 


Although the first term above is singular at a; = — cos it becomes finite when combined with the term proportional 
to Fc (^y/y(J + a;-)^ To see this note that 
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sin^Vl — ■ 


|m| 


r^iTTl4: 


I cos c + 2^1 \ 1 + X cos C + |x + cos ^1 I 
1 I sin CVl — 


+ u 


TT sin CVl — x^ v^l + X- 

\ \m\ 


F 


(v^Ri 


1 


cosC + x| \^1+xcosC+|x + cosC|y ^2 sin Cv^F^Vl + x_ 

\ Nl 


+ x_ 


+' 


(G13) 


1 


sinC\/l — : 


1 + x+ 


cos C + x| I I 1 + X cos C + |x + cos Cl / Y 2 sin CVl - x^ 


where we used 


+ X+ ^J\ + X- = |x + cos Cl , 


(G14) 


to get the last line, and where the dots correspond to the Fresnel cosine and sine integral terms from F. Since, to 
leading order in x + cosC, the expression in curly brackets is —|m|| cosC + x|/sin^ C) it follows that (G12) for Im{y,x) is 
actually finite at x = — co s C and therefore integrable. For small values of the argument C{y) « y and S{y) « 7ry^/6, 
so the terms in Eq. (G13l represented by the dots are also finite for all x, and proportional to ^/y near x = — cosC- 
In deriving expression (G12), we have neglected some terms of 0(l/y^), but terms of that order and higher are 
present in Eq. (G12) so these orders have been treated inconsistently. To obtain a consistent result at 0(2/i, 1 / 2)1 
could expand this expression and drop terms of higher order. However, keeping the incomplete higher order terms in 
Eq. (G12) was found empirically to give a better approximation to numerically computed overlap reduction functions. 


Appendix H: Analytic calculation of the overlap reduction function for co-directional pulsars for 

scalar-longitudinal backgrounds 

For two pulsars that lie along the same line of sight as seen from Earth (i.e., cosC = 1), the calculation of Im{y,x) 
can be done analytically. For this case 


Im{y,x) 
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COS(^ = l Jq 
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The integral for then takes the form 




/-I 

tL/ 


(1 -\-x)^ 

= nNrSmo [Gfiyi) + Gfi-y^) - Cf (yi - ^2)] , 


where 


Gf{y)= f^dx 


By making the expansion 


(1 + xy 


= 3 — 2 x + x^ — 


(1+x) {l+xY 


(H2) 


(H3) 


(H4) 


we can write 
Gf{y)= j\x 


Z — 2 x -\- x^ — 


(1+x) {l + xY 

2{-iye 

- dHi{y) + Ki{y) 


20 r 4 ^ 4 , -iy 

= y <5zo - -^Sii + —5/2 - 2(-i)'e ^ 


Pi{x) (1 - 


“ 3 ) ji{y) + 


y" 


y 


21 + 1 

y 


+ 2i] ji+i{y) - ji+ 2 {y) 


(H5) 


where 


Hi{y) = J (1 ® * 2 /( 1 +“)^ , 

Ki{y) = da/ (l - . 


(H6) 

(H7) 


These last two integrals are most easily evaluated using the recursion relation in Eq. (D7) for Legendre polynomials, 
which for this calculation is most conveniently rewritten as: 


Piix) = - ^—j^Pi- 2 ix) + ^—-j^{l+x)Pi-i{x), ioT l> 2. 


(H8) 


This leads to 


Ho{y) = Cin(2y) + iSi(2y), 

H^{y) = -Ho{y) + 2 +^-{l-e-^^y) , 

Hi{y) = -^?^i7/-i(y) - ^^i7/_2(y) - 2(-a)'-i^^^e-*^j/_i(y), for Z > 2, 


Ko{y) = ^ +ySi(2y)^ +i sin(2y) - y Cin(2y) + y 

Ki{y) = Ho{y) - Ko{y ), 

Ki{y) = -^^^K,-,{y) - i^Ki.2{y) + ^?^i7/-i(y), 

where Si(a;) and Cin(a:) denote the sine and cosine integrals respectively, defined by 

1 — cos t 

dt . Cin(a;) = j dt 
/o 


1 + 


da; 


Li 1 + a;_ 
for Z > 2 , 


Si(a;) = f dt —— . Cin(a;) = f 
Jo t Jo 
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(H9) 


(HIO) 


(Hll) 


Note that the last two terms (in square brackets) in the above expression for K[){y) will cancel when forming the 
combination iLo(yi) + Koi~y 2 ) — Koiyi ~ 2 / 2 ); which enters the expression for r^(/). The above recursion relations 
are particularly useful when the values of 77/(y) and 77/(y) are required at fixed y for all I < /max- 
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For isotropic backgrounds (/ = m = 0), an expression for the scalar-longitudinal overlap reduction function for 
equidistant (yi = 1/2 = y)-, co-directional (cos^ = 1) pulsars valid in the limit y ^ 1 was given in Chamberlin and 
Siemens |34| . Equation (H2) reduces to that result in the appr opri ate limit, as we now show. 

For equidistant pulsars and I = 0,m = 0, the last term in Eq. (H2 1 is Gq (0), which is zero. This can be seen by direct 
evaluation or by noting that the last term in square brackets in Eq. (H51 reduces to [3jo( 2 /) + (l/j/ + ‘^i)ji{y) — j 2 (y)] 
for 1 = 0, which tends to 10/3 as y —?► 0. When multiplied by the pre-factor of —2, this cancels the first term in 
Eq. (H5). Likewise, i?o(0) = 0 and Kq{0) = 0, so Gg (0) = 0. The equidistant, co-aligned, isotropic overlap reduction 
function is therefore 


pi 

^ 00 


(/)|cosC=i — ^ [Gg (y) + Gg (-y)] — — [Gg (y) -I- Gg (y)*] — ^/tt Re{Gg (y)} . 


(H12) 


We now evaluate this expression in the limit y 3> 1. All spherical Bessel functions decay to zero as y —>■ 00 , so the 
term in square brackets in Eq. (H5| makes no contribution in this limit. Hence, we focus on the behaviour of i?o(y) 


and Aro(y) for large y. We make use of the following asymptotic expressions: 

Si(y) « I, y > 1 , 

Cin(y) « 7 -f ln(y), y » 1, 

where 7 is the Euler-Masheroni constant, 7 = 0.57722 • • •. We deduce that, for large y, 


20 
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Gg (y) « ^ - 4(7 -f ln( 2 y)) -i2Tr- ^ + ^ - iy{j + ln( 2 y)) -h ^e 
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(H13) 


(H14) 


SO 


and 


37 TTU 1 

Re{Go (y)} « y - 47 - 41 n( 2 y) + y + 2 cos( 2 y), 

« y-47-41n(2y)-hy , 


(H15) 


rgg(/)|cosC=l ~ 


y -47-41n(2y) -f y 


37 , /47r/L 

-- -47-41n - 

6 \ c 


TT^fL 


(H16) 


where / is the gravitational-wave frequency and L is the distance of the two pulsars from the Earth. This agrees 
with Eq. (40) of Chamberlin and Siemens [34], apart from a factor of y/n, which comes from a difference in our 
normalization convention. 


Appendix I: Analytic calculation of the overlap reduction function for anti-directional pulsars for 

scalar-longitudinal backgrounds 


For two pulsars that lie in antipodal positions along the same line of sight as seen from Earth (i.e., cosC = —1), 
the calculation of Im{y,x) can also be done analytically. For this case 


cosC=-l 7o 1 — xV / l — x\ J 


The integral for r^^(/) then takes the form 
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we can write 

rL(/) 


cos —1 


= TTN^5mO / da; 


/-I 
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Pi 


= nNrSmO 


-,,Sio--Si2 + 2{-tye-^y^ 
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yi 


y2 
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2 / +1 

ji(yi + y 2 ) H-^— ji+iivi + y 2 ) - ji+ 2 (yi + 2 / 2 ) 


2/1 + 2 / 2 ' 


+^Hi{yi,y2) + ^Ht{y2,yi) 


where 


Hi(yi,y2) = y'^d;r^^^P;(a;)(l-e-**'di+-)) _ 

This final integral can be obtained via a recurrence relation using Eq. ( |H 8 [ ) from App. We find 
Ho(yi,y 2 ) = Cin(2yi) + iSi(2yi) + (Cin(2yi) + 2Si(2j/i) - Cin [2(yi + j/ 2 )] - *Si [2(yi + j/ 2 )]) 

sMyi + y2) 

yi + 2/2 


Hi(yi,y 2 ) = -^ 0 ( 221 , 2 / 2 ) + 2 ( 1 - 


2/1 


y2 


Hi(yi,y 2 ) = - ^^—^Hi-i(yi,y 2 ) - gi- 2 ( 2 /i, 2 / 2 ) 




e *^^ji-i( 2 /i)+e*^=jz_i(?/ 2 ) - v^)ji_^(yj^ + y^)] , for 2 > 2 , 


(14) 


(15) 


(16) 


where, as before, Si(x) and Cin(a;) denote the sine and cosine integrals, which were defined in Eq. (Hll). The result 
for Ho(yi,y 2 ) can be obtained by rewriting Eq. (|l5|) as a combination of integrals of the following four forms: 


r 2 y 
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('■2y 


du 


1 — cos u 
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1 — cos u 
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cos(aii) = icin[2(a + 1)?/] + icin[2(a — l)y\ — Cin(2ay), 
sin(au) = Si(2ay) - ^Si[2(a + l)y] - ^Si[2(a - l)y], 


r 2 y 


r 2 y 


Qi-n 7/ 1 1 

du -cos(aw) = -Si[2(a + l)y] — -Si[2(a — l)y] , 

u 2 2 

sin 7/1 1 

du -sin(aw) = -Cin[2(a + l)y] -Cin[2(a — l)y]. 

u 2 2 


(17) 


Appendix J: Analytic calculation of the overlap reduction functions for vector-longitudinal backgrounds 


Ignoring the pulsar terms, the overlap reduction functions for an uncorrelated, unpolarised, anisotropic vector- 
longitudinal background are given by I'Ymif) — d and 
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where 


Sin 

^Im — —l) / dx X 
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-1 '■ ■' 1 + X 
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1 + X cos C + a/ 1 — x^ sin ( cos (j) 


The integral Kim{x) can be evaluated using contour integration, as described in for the response of a PTA to 
anisotropic backgrounds with GR polarisations. The result is 
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The non-zero Iim’s can be straightforwardly evaluated: 
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The Jim’s can be written in terms of the m(cosC) functions defined in 


r,L,m(cos C) = 


Ptr,L,micOsC) = 
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l-cosC i^ + xYdx-- 
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while for to > 0 we have 
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LH 
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and = (—1)”’'A^/'" Jzm- Explicit expressions for the J'^rLm(cosC) functions are given in App.j^ 


(J7) 


1. Limiting case: cosC = 1 

As noted in the main text, in the limit cosC —t 1, the to = 0 overlap reduction functions calculated above diverge. 
This singularity is eliminated if the pulsar terms are included in the integrand, and the pulsars are assumed to 
be at finite distance. Proceeding in a fashion identical to the case of co-directional pulsars in scalar-longitudinal 
backgrounds, we find 
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where 


G^{y) = dx (l - e-*^(^+")) . 




—2 -|- 2 ic — X “h 
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14 4 4 ,- 
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^+ 2 . 1 - 2 ),(.) 
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(J9) 


+ ‘2Hi{y), 


with Hi{y) defined as in Eq. (H6). This is a finite expression provided yi and y 2 are finite. 


Appendix K: Evaluating the integrals for transverse tensor and vector-longitudinal backgrounds 

The overlap reduction functions for both the standard transverse tensor and vector-longitudinal backgrounds can 
be written in terms of the functions 
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These integrals can be evaluated using the series representation of the Legendre polynomials 
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for which 
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For the standard transverse tensor backgrounds, we also need to evaluate -7^ j. i C) for r = — 1. 
reduced to combinations of /"“g ; ^(cos^) and g i writing (1 — cc) = 2 — (1 + x): 

= 2J'-g ; ,„(C0SC) - -^■+l,o,i,m(cOsC) ■ 


This can be 


(K8) 


Alternatively, we can just evaluate this integral directly, finding 
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Appendix L: Recovering the overlap reduction function for an uncorrelated, anisotropic scalar-transverse 

background 

Ignoring the pulsar term, we can show that the response of a pulsar to the indiviudal modes of a scalar-transverse 
gravitational-wave background can be use d to recover the overlap reduction function for an arbitrary uncorrelated, 
anisotropic background. Inverting Eq. (621 to find aE)(/) gives 


4m)U) = V2 / d^fl^ hB{f,k)Y:^{k) , 


IS2 


(LI) 


from which we deduce the following quadratic expectation values: 
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where C'imi'm'(/>/O = ^Eni'm'^B{f)5[f — f) assuming stationarity. For a Gaussian-stationary, uncorrelated, 
anisotropic background, the quadratic expectation value of breathing mode amplitudes is given by Eq. (43|: 


{HbU, k)h*B{f, k')) = HB{f)PBCk)Syk, k')S{f - f). 


(L3) 


The angular distribution of gravitational-wave power can be expanded in terms of scalar spherical harmonics (see 
Eq. (44)). Hence the integrals over the sphere in Eq. (L2) can be explicitly evaluated: 
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Now the overlap reduction function between pulsars 1 and 2 is given by 


p-B _ \ ^ \ ^ /^B tdB j:>B^ 

^ / V / V ^lrnl'rn'^l{lrri)^2{l'rn’) 

{Im) {I'm') 
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jfB pB 
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0 0 0 ) Pl(lm)^2(l'm') (L5) 


Note that the breathing response is limited to Z = 0,1. Hence, Wigner-3j selection rules restrict the sensitivity 
of the bre ath ing mode over lap reduction function to L < 2. By substituting the breathing response fu nction 
from Eq. (|6l|) into Eq. (L5), we fully recover the form of obtained by direct calculation in Eq. (481. For 

example, with L = 0,M = 0, Eq. (L5) gives E^q = (•y/7r/2)(l + | cosC), where ^ is the angular separation be¬ 
tween the two p ulsars. (Recall that = •\/4^/2 for an isotropic uncorrelated background, as described at the 
end of Sec. IVE ) This exactly matches the expression given by Eq. (481, as do the remaining expressions for L = 1, 2. 
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